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1 Introduction 



The masses of the light quarks are not directly accessible to experiment and have to 
be determined through a non-perturbative calculation, taking low-energy hadronic data 
(such as the kaon masses) as input. Chiral perturbation theory is able to predict their 
ratios with a precision of a few percent M. In order to obtain estimates of individual 
quark masses, it thus suffices to determine a particular linear combination using lattice 
QCD g-|l| or QCD sum rules @-[2§. 

Calculations of quark masses in lattice QCD are in principle straightforward, but 
there are several sources of systematic errors which must be studied carefully (for re- 
cent reviews of the subject sec [24,25]). One of the principal uncertainties arises from 
the renormalization constant needed to convert from lattice normalizations to the MS 
scheme of dimensional regularization. Usually this factor is only known in one-loop 
perturbation theory in the bare coupling. The limitations of lattice perturbation theory 
are well known, and in order to remove all doubts about the reliability of quark mass 
calculations in lattice QCD, it is evident that a non-perturbative determination of the 
renormalization factor is required. 

The general problem of non-perturbative renormalization of QCD addresses the 
question how the perturbative regime of QCD is related to the observed hadrons and 
their interactions. This relation involves large scale differences, and thus the problem 
is not easily approached using numerical simulations, for which only a limited range of 
lattice spacings is accessible. 

A strategy how to overcome this technical difficulty has been described in [26], and a 
comprehensive introduction into the subject is presented in [27]. An alternative method 
was introduced in |2£| and is actively being pursued [|, 12, 29, 30]. In our strategy 
the key idea is the introduction of an intermediate renormalization scheme, with a 
controlled (non-perturbative) relation to the lattice normalizations, and in which the 
scale evolution of masses and couplings can be computed non-perturbatively from low 
to very high energies. In the high-energy regime one can continue the scale evolution 
using perturbation theory, which allows for the determination of the renormalization 
group invariant quark masses M and the A-parameter. Since the matching of M and A 
between different renormalization schemes is exactly computable, it is evident that all 
reference to the intermediate scheme is eliminated at this stage. 

Thus, the problem of non-perturbative quark mass renormalization can be split 
into two parts. The first is the computation of the scale dependence of the running 
masses fn in the intermediate scheme from low to high energies and its relation to the 
renormalization group invariant quark masses M. The ratio fn/M is then known at 
all energy scales covered by the calculation. The second part is the matching of the 
renormalized masses fn with the bare current quark masses m(go). This amounts to the 
computation of fn/m(go) at a certain value of the energy scale which is chosen as the 
matching point. Since the matching can be performed at low energies, no excessively 
large lattices are required to contain all relevant scales. Whereas the second part clearly 
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depends on the details of the lattice regularization, such as the chosen lattice action and 
the bare coupling, the ratio fn/M is a universal quantity. 

In this work we have calculated the regularization independent factor fn/M in 
quenched QCD, using the Schrodinger functional (SF) as our intermediate renormaliza- 
tion scheme. We emphasize that the SF scheme is chosen entirely for practical purposes, 
since it allows the non-perturbative computation of the scale dependence over several 
orders of magnitude in a controlled way, using a recursive finite-size scaling technique 
and numerical simulations. The numerical results for the scale dependence of fh can be 
extrapolated reliably to the continuum limit, so that a truly universal result for fn/M 
is obtained. 

In addition, we have computed the matching between the renormalized quark 



masses in the SF scheme and the bare quark masses in 0(a) improved lattice QCD [31] 
for a range of bare couplings at a fixed scale. Thus, by combining the universal ratio 
fn/M with the matching factor fn/m(go) we have obtained the total renormalization 
factor M/m(go), which relates the bare current quark masses to the renormalization 
group invariant quark masses. 

We also report on our determination of the A-parameter in quenched QCD. For this 
purpose, we have supplemented the numerical data for the running coupling published in 
|32j |, so that the scale evolution could be traced over more than two orders of magnitude. 
Preliminary results of our calculations have been presented in p3[ . 

The rest of the paper is organized as follows. In Sect. ^ we discuss the scale depen- 
dence of the running coupling and quark masses in perturbation theory. The Schrodinger 
functional scheme is described in detail in Sect. ||. The step scaling functions, which de- 
scribe the scale evolution of the running parameters, are briefly discussed in Sect. ||. In 
Sect. H we describe the extraction of fn/M and the A-parameter. The matching of these 
results to hadronic schemes is presented in Sect. [| and Sect, ^contains our conclusions. 
In order to avoid distracting the reader from the main results, we defer all technical 
details about the computation of the step scaling functions, the error propagation in the 
scale evolution and the calculation of fn/m(go) to Appendices [A|, |B| and |C], respectively. 

2 Running coupling and quark masses in perturbation theory 
2.1 Mass-independent renormalization schemes 

QCD is a theory with iVf + 1 free parameters: its coupling constant g and the masses 
of the JVf different quark flavours, {m s , s = 1, . . . , N{}. In the definition of the theory 
through a regularization, these are initially taken to be the bare parameters in the 
Lagrangian. Upon removal of the regularization, one defines renormalized parameters, 
g, {fn s ,s = 1, ...,iVf} at some energy scale [i. In the following we assume that the 
normalization conditions which are imposed to define g, fn s are independent of the quark 
masses themselves. Examples for such mass-independent renormalization schemes are 
the MS scheme of dimensional regularization [34, 33] and the Schrodinger functional 
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(SF) scheme [26, 36|. The latter will be described in Sect. ||. The renormalized running 
parameters are then functions of the renormalization scale [i alone, and their scale 
evolution is given by the renormalization group equations 

dg_ 

H^r-^- = T(q)m x , S = 1, 



r(g)rn s , 



(2.1) 
(2.2) 



The renormalization group functions, (3 and r, have a perturbative expansion 

.'< J 7. , 7, ..2 | 



Kg) s - 



-g 



[b + hg 2 + b 2 g 4 + 



r(g) S ~° -g 2 {d Q + d 1 f + d 2 f + ..) . 
The universal coefficients bo, b\ and do are given by 

b = (4vr)- 2 {ll - fiVf} , &! = (47T)" 4 {102 - f Nt} , 
d = 8/(4tt) 2 , 



(2.3) 
(2.4) 

(2.5) 
(2.6) 



whereas the higher order coefficients b 2 , b^,... and d±, d 2 ,... depend on the chosen 



renormalization scheme. Equations (2^,2^) can be solved straightforwardly. This 
yields the exact relations between the scale dependent quantities g(/j,),m s (/j,) and the 
renormalization group invariants, 



A 



(bof r^e-- 1 ^ 2 x 



cxp 



dg 



1 1 

+ 



P(g) b g 3 b 2 g 



(2.7) 



M, 



m s (2b g 2 r do/2bo x 



exp 



dg 



r(g) _ do_ 
0(g) b g 



(2. 



Unlike the case of the A-parameter, there is no universally accepted normalization con- 
vention for the renormalization group invariant quark masses M s . Eq. (2.8) complies 
with the conventions used by Gasser and Leutwyler . 

In contrast to g and m s , the renormalization group invariants A and M s do not 
depend on the scale; from a technical point of view they are the integration constants 
of the differential equations (2.1,2~2|). From eq. ( [2. &] ) and the scale independence of 
M s , we read off immediately that ratios of running quark masses, fn s /fn s i, are scale 
independent in mass-independent renormalization schemes, and are equal to M s /M s i. 

The renormalization group invariant quark masses are easily shown to be inde- 
pendent of the renormalization scheme, while the A-parameter changes from scheme to 
scheme by factors that can be calculated exactly. It is hence natural to take A and 
the renormalization group invariant quark masses M s as the fundamental parameters 
of QCD. 

Furthermore one concludes that in order to extract the fundamental parameters of 
QCD, one may choose a scheme which is particularly suited for the computation of A 
and M s . 
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Figure 1: Running coupling a-^jg = <?|jg/47r in quenched QCD. From the dotted to the full curve, the 
perturbative accuracy increases from two to four loops. 



2.2 Illustration: the MS scheme 



In the MS scheme the perturbative expansions eqs. (2~lj|,2,4) are known to four loops 
lH-@. Here we use them to plot the evolution of the running coupling and quark 
masses in this scheme for iVf = flavours. The dotted, dashed and full curves in Figs. |l] 
and |2| have been obtained from eqs. (2/7) and ( f?lj| ) by substituting the 2-, 3- and 4- loop 
expressions for the (3- and T-functions and computing the integrals in the exponents 
numerically. 

Conventionally the running coupling and quark masses are quoted in the MS scheme 
at a certain value of the normalization mass fi. Figures [l] and [2] illustrate that, once A 
and M are known, the running parameters in the MS scheme are uniquely determined 
to any given order of perturbation theory. 

As we shall see in Sect. |6| the value of A^jg in quenched QCD is about 240 MeV. 
From Figs, [l] and [2] one infers that perturbation theory appears to "converge" at energies 
as low as 1 GeV. However, since the MS scheme is only defined to any finite order of 
perturbation theory, it is impossible to make a solid statement about the total uncer- 
tainty in the running coupling and quark masses at a certain reference scale. In other 
words, the running coupling and quark masses in the MS scheme are only meaningful 
to a given order of perturbation theory. 
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Figure 2: Running quark mass in quenched QCD. The flavour index has been omitted, since the entire 
graph is independent of which quark flavour is considered. 



3 The Schrodinger functional scheme 

In this section we describe the intermediate renormalization scheme, which we have 
used to relate hadronic observables to the renormalization group invariant parameters. 
The chosen scheme, defined using the Schrodinger functional, has been the subject of 



a series of publications [p3, 36, For clarity we repeat the basic concepts of the 



Schrodinger functional but shall otherwise be rather brief and refer the reader to the 
literature for details. A pedagogical introduction can be found in 27]. 

3.1 General definitions 

As already mentioned in the introduction, the main difficulty of relating the observed 
hadron spectrum to the perturbative regime of QCD in a controlled way is the large scale 



difference. Non-perturbative renormalization such as the approach introduced in [28] 
relies on all relevant physical scales being accommodated on a single lattice. This results 
in a rather narrow energy range which can be explored. 

It has been demonstrated [^J that this problem can be overcome by simulating 
a sequence of lattices with decreasing lattice spacing. Any single lattice only covers a 
limited range of distances, but through the use of a finite- volume renormalization scheme 
it is possible to match subsequent lattices. With only a few steps of this procedure one 
can easily cover a much larger energy range than otherwise possible. 
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The Schrodinger functional is a particular finite- volume renormalization scheme. It 
is based on the formulation of QCD in a finite space-time volume of size Txl 3 , with 
periodic spatial boundary conditions and inhomogeneous Dirichlet boundary conditions 
at time xq = and xq = T 49|, 5(|. This is realized by requiring the spatial components 



of the gauge field at the boundary to be equal to some prescribed constant abelian fields 
C and C . This choice of boundary condition introduces a frequency gap on the quark 
and gluon fields, so that simulations for vanishing quark mass can be performed. 

We take over the exact form of the action and the boundary conditions from 



Sect. 2.1 of ref. [52]; any notation not explained here is taken over from that refer- 
ence. The definition of the scheme requires that the ratio T/L is kept fixed. For reasons 
that will be explained below we have set T = L throughout this work. Within this 
set-up, renormalization conditions are then specified at scale fj, = 1/L and vanishing 
quark mass. Thus, the Schrodinger functional is a mass-independent renormalization 
scheme, in which the running coupling and masses scale with the box size L. 

3.2 Running coupling and masses 

Our next task is to specify exactly the definition of the running coupling and quark 
masses in the SF scheme, which we have used in the numerical simulations. In particular, 
we will exploit the freedom in choosing the boundary fields C and C , and the angle 9 
which appears in the spatial boundary conditions of quark fields. From a practical point 
of view, these have to be chosen with care in order to ensure that the scale evolution 
of masses and couplings can be determined with high accuracy in the continuum limit 
and that contact with the perturbative evolution at high energies can easily be made. 
The "optimal" definitions of g and rh [32, 36, 4£, 52] are obtained by considering the 



following criteria: 

• Choose observables with small variance in Monte Carlo simulations to ensure good 
statistical precision of the results; 

• Avoid parameter values which introduce "accidentally" large coefficients in the 
perturbative expansion of e.g. the (3- and T-functions; 

• Minimize the contamination of coupling and quark mass by lattice artefacts. 
The renormalized coupling at length scale L is given in terms of an infinitesimal 



variation of the boundary conditions for the gauge fields, exactly as defined in refs. [32, 
[52| , where all details can be found. It is denoted by g(L), whereas in accordance 
with standard notation in the MS scheme, we choose the corresponding energy scale as 
argument for a, viz. 

a(») = ^-, » = 1/L. (3.1) 
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We now describe the definition of the running quark mass m(/z) in more detail and 
introduce it formally in the framework of continuum QCD. Its definition is easily made 



rigorous in the lattice regularization (see Appendix A. 2 and pq| ). The starting point is 
the PCAC relation, 



(3.2) 



(3.3) 



(3.4) 



Once eq. (|3.2| ) is chosen to define the quark masses, their normalization is given in terms 
of the normalization of A R and -Pr. The axial current is normalized naturally through 
current algebra relations 54-S2| and does not acquire a scale dependence through renor- 
malization. Therefore, both scheme- and scale-dependence are contained in the defini- 
tion of Pr. 

In order to formulate a normalization condition for the pseudoscalar density, we 



d^A^fj, = (m s + m s ,)P R , 
between the renormalized axial current, 

and the pseudoscalar density, 

P R (x) = Z P i(j s (x)j 5 ip s/ (x) . 



set 1 36] 



L, C = C' = 0, 6 = 1/2. 



(3.5) 



This choice is motivated by the criteria listed above, and a more detailed explanation 
is given in ref. [p6|. We then start from the correlation functions 



fp(x ) = -iy dVd d z(V(a:)75^V(x)C(y)75^ a C(z)) (3.6) 
fx = -^6/d 3 ud 3 vd 3 yd 3 z(C , (u) 75 ir a C / (v)C(y)75^C(z)), (3.7) 



involving the boundary quark fields, Cj • • • > s [63]- Here and in the following, the Pauli- 
matrices, r a , are understood to act on the first two flavour components of the quark 
fields. As the renormalization is flavour-blind in a mass-independent renormalization 
scheme, this is sufficient to define the renormalization constants. Furthermore we in- 
troduce bare current quark masses m s via 



^s{ x )l)J,lh^s'{x) 



and define the renormalization constant 



Z P (L) 



fp(L/2) ' 



at 777, 



0, 



(3.8) 



(3.9) 
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Here, the numerator yffi is introduced to cancel the multiplicative renormalization of 
the boundary quark fields £ and which appear in the definition of /p. Our convention 
eq. (3.9) ensures that Zp = 1 at tree level of perturbation theory. 



Explicitly the definition of the running quark masses in the SF scheme now reads 

Z\ 

m(fJ,) s = Zp ( L ) m s> V=l/L. (3.10) 
The above expressions are easily given a precise meaning in the lattice regulariza- 



tion, and we refer to the formulae listed in Appendix A. 2.1. We emphasize, however, 
that the SF-scheme is not linked to a particular regularization. Note that Zp is defined 
in terms of correlation functions at finite physical distances for which on-shell 0(a) 
improvement may be applied [63, 64] to reduce lattice artefacts. Although this detail 
is irrelevant for the definition of the scheme itself, it will be important in its numerical 
implementation. 

From now on we will drop the flavour index s on all quark masses, except when we 
refer to a particular flavour. 

Finally we note that the normalization constant Za of the axial current is known 
non-perturbatively in 0(a) improved quenched lattice QCD |52]. Hence, in order to 
determine the complete renormalization factor in eq. ( |3.10| ), it suffices to compute the 
renormalization constant Zp. 

3.3 Perturbation theory 

The SF scheme has been studied perturbatively in quite some detail. In the rest of 
this section we list the results which are relevant for the present computation, setting 
Nf = 0. All quantities where we do not indicate the scheme explicitly refer to the SF 
scheme. The A-parameter is translated to the MS scheme via |32]| 

A = 0.48811(l)A M g , (3.11) 



and the three- loop coefficient of the (3- function was determined as [65| 

b 2 = 0.483(9)/(4vr) 3 . (3.12) 
In addition, the two-loop anomalous dimension of the quark mass, 



di = 0.217(l)/(4vr) 2 , (3.13) 



was computed as part of this project |36|. Given the perturbative expansion coefficients 
of (3 and r in the MS scheme, these equations fix the one-loop relation of fn to ??%jg 
and the two- loop relation of a to ct™. Explicit formulae may be found in [36, 65fl. 



For our purposes, the coefficients 62 and d\ serve primarily to continue the scale 
evolution of the running coupling and quark mass in the SF scheme to infinite energy 
using perturbation theory. This allows for the extraction of A and M, which can be 



easily converted into their counterparts in the MS scheme through eq. ( 3.11 ) and the 



fact that M = M^g. 
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u 


a(u) 


a P (u) 


0.8873 


0.981(9) 


0.9683(21) 


0.9944 


1.110(11) 


0.9672(23) 


1.0989 


1.252(11) 


0.9622(24) 


1.2430 


1.416(16) 


0.9579(29) 


1.3293 


1.528(20) 


0.9470(28) 


1.4300 


1.703(24) 


0.9407(30) 


1.5553 


1.865(23) 


0.9382(33) 


1.6950 


2.095(25) 


0.9297(32) 


1.8811 


2.412(32) 


0.9284(36) 


2.1000 


2.771(41) 


0.9168(37) 


2.4484 


3.464(40) 


0.8942(38) 


2.7700 




0.8781(42) 


3.4800 




0.8451(55) 



Table 1: Simulation results for the step scaling functions a(u) and crp(ii) 

4 Step scaling functions 

In the SF scheme a change of the renormalization scale amounts to a change of the box 
size L at fixed bare parameters. By considering a sequence of pairs of volumes with 
sizes L and 2L, one can thus study the evolution of the running coupling and quark 
masses under repeated changes of the scale by a factor 2. Effectively one constructs a 
non-perturbative renormalization group in this way. 

The evolution from size L to 2L of the running coupling and the normalization 
of the pseudoscalar density is described by the step scaling functions a(u) and ap(u) 
according to 

g\2L) = a(u), u^g 2 (L), (4.1) 
Z P (2L) = ap(u)Z P (L). (4.2) 

For a given lattice resolution a/L both functions can be computed non-perturbatively 
through numerical simulations of the Schrodinger functional. It is important to realize 
that the box size L can be as small as we like in these calculations. The only restriction 
is that L and the low-energy physical scales in the theory should be significantly larger 
than the lattice spacing to avoid uncontrolled cutoff effects. 

The step scaling function associated with the coupling has first been computed 
in ref. |3^]. We have extended these calculations so that u(u) is now known precisely 
for altogether 11 values of u (see Table [l]). Moreover we have obtained ap(u) at the 
same couplings and, in addition, at u = 2.77 and u = 3.48. The technical details 
of our computations are summarized in Appendix |A[ Here we only mention that the 
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Figure 3: Polynomial fit of the data for the step scaling function a(u). The errors on the data are 
about equal to the symbol size. 

lattice results for the step scaling functions depend on the lattice spacing and have to 
be extrapolated to the continuum limit. The extrapolation is required but not critical, 
since in all cases considered the observed lattice effects are small, particularly after 
including the relevant 0(a) correction terms. 

5 Running coupling and quark mass in the SF scheme 

Having computed the step scaling functions, we can now work out the evolution of 
the coupling and the quark masses over roughly two decades of the energy scale. It 
should be stressed from the beginning that all results obtained in this section refer to 
the continuum theory. The basic idea is to start from some initial values at low energies 
and to step up the energy scale by factors of 2 by repeated application of the step scaling 
functions. In terms of the box size L this leads us from some large size L max to smaller 
sizes 2~ k L max . Eventually the perturbative regime is reached and perturbation theory 
may then be used to extract the A-parameter and the renormalization group invariant 
quark masses. 

The largest value of g 2 which can be reached with the available data for the step 
scaling function a(u) is 3.48. We thus define the scale L max through 

g\L)=3A8 at L = L max . (5.1) 
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Figure 4: Comparison of the numerically computed values of the running coup ling in the SF scheme 
with perturbation theory. The dotted and dashed curves are obtained from eq. (2.7) using the 2- and 
3-loop expressions for the /3-function. The errors on the data are smaller than the symbol size. 



The sequence of couplings^ 

u k = g 2 (2- k L max ), k = 0,1,2, ... , (5.2) 

can then be calculated by solving the recursion 

uq = 3.48, (j(u k+ i) = u k . (5.3) 

A technical problem here is that cr(u) is only known at certain values of the coupling and 
only to a finite numerical precision. This difficulty can be resolved by fitting the data 
with a polynomial, as shown in Fig. |||, and using the fit function in the recursion (5.3). 



The analysis presented in Appendix |B] shows that the systematic uncertainty which is 
incurred by this procedure is negligible compared to the statistical errors, provided one 
stays in the range of couplings covered by the data. 
After six steps the recursion yields 

f{L) = 1.053(12) at L = 2~ 6 L max (5.4) 

and two more iterations bring us down to a value of 0.865(11). As will become clear 
below, these couplings are deep in the high-energy regime where the scale evolution is 



2 Note that, contrary to the convention used in (33], the box size L decreases for increasing k. 
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Figure 5: Polynomial fit of the data for the step scaling function op (it). 



accurately described by perturbation theory. The A-parameter may thus be calculated 
by setting fj, = 1/L in eq. ( pT7| ) and inserting the coupling quoted above. To evaluate 
the integral in this formula, one may safely use the 3-loop expression for the /3-function 
since the coupling is small in the whole integration range (if an approximately geometric 
progression of the coefficients bk is assumed, the error from the neglected higher-order 
terms is an order of magnitude smaller than the statistical error) . As result one obtains 

A = 0.211(16)/L max (5.5) 

and from Fig. ||] one can now see that the numerically determined values of the coupling 
are indeed closely matched by the perturbative evolution in the range where the A- 
parameter has been extracted. In this plot the error in eq. ( |5.5| ) has not been taken 
into account since it amounts to an overall scale shift which is of no relevance for the 
comparison of the scale evolution of the coupling with perturbation theory. 

We now proceed to discuss the scale dependence of the renormalized pseudoscalar 
density and the running quark mass. In this case the data for the step scaling functions 
listed in Table [l] allow us to trace the scale evolution up to L = 2L max . If we define 

v k = Zp(2- fe+1 L max )/Z P (2L max ), k = 0, 1, 2, . . . , (5.6) 

the recursion to be solved is 

vq = 1, v k +i = v k /ap(u k ), (5.7) 
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Figure 6: Comparison of the numerically computed values of the running quark mass in the SF scheme 



with perturbation theory. The dotted, dashed and solid curves are obtained from eqs. (2.7) and (2 
using the 1/2- , 2/2- and 2/3- loop expressions for the r- and /3-functions respectively. 



where u^, k > 0, is the sequence of couplings introduced above. Using the fit for ap(u) 
displayed in Fig. [5] it is straightforward to perform this calculation and after seven steps 
one finds 

Zp(L)/Z P (2L max ) = 1.759(15) at L = 2~ 6 L max (5.8) 

(see Appendix || for further details). At this point the coupling is so small that it is 
safe to insert the two- and three-loop expressions for the r- and /3-function in eq. fl2.8|) . 
Taking eq. (|5.8|) into account the result 



M/m(p) = 1.157(12) at 11 = (2L n 



(5.9) 



is then obtained. Similar to the coupling, the scale evolution of the running quark masses 
in the SF scheme is accurately reproduced by perturbation theory down to surprisingly 
low energies (see Fig. ||). 



6 Matching to hadronic observables 

At this point the A-parameter and the renormalization group invariant quark masses M 
are known in terms of the reference scale L max and the running quark masses fn(n) at 
[i = (2L max ) _1 . To complete the calculation we now need to relate L max and m(/x) 
to a hadronic scheme, where the parameters of the theory are fixed by requiring a set 
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of low-energy hadronic quantities to assume their physical values. It should again be 
emphasized that this step does not involve any large scale differences. In particular, 
the calculations can be carried out using lattices which cover all relevant scales, as in a 
conventional hadron mass computation. 

6.1 Computation of 

We now describe how to convert the result for the A-parameter in eq. ( |5.5| ) into A^g 
expressed in terms of a hadronic scale. A convenient hadronic reference scale is the 
radius ro which has been introduced in ref. (6(|. It is defined through the force between 
static quarks and is thus a purely gluonic quantity in the quenched approximation, 
similarly to g 2 (L) and L ma x- Relating L max to ro is straightforward, but requires data 
for ro/a with good precision for a range of go, i n order to obtain the continuum limit [67]. 
As part of our project, the ratio L m axAo was computed 



W/r = 0.718(16) . (6.1) 
We can now remove all reference to the intermediate renormalization scheme and quote 



the A-parameter in units of ro- By combining eqs. fl6.1|) , (|5.5[) and ( 3.11 ) we obtain 



A<| = 0.602(48)/r , (6.2) 

where the superscript reminds us that this quantity has been determined for Nf = 0, 
i.e. in the quenched theory. For illustration, physical units can be introduced by setting 
ro = 0.5 fm, which gives 

= 238(19) MeV. (6.3) 

However, as is well known, the conversion to physical units is ambiguous in the quenched 
approximation. The solid result is eq. (|6.2| ). 

6.2 Computation of the total renormalization factor 

Next we describe the calculation of the renormalization factor Zyi(go) which relates the 
bare current quark masses m to the renormalization group invariant quark masses M. 



The renormalized current quark masses have been defined in eq. (3.8), and their def- 



inition in the 0(a) improved lattice theory is given in eq. (A.7). By including the 
appropriate renormalization factors for the axial current and the pseudoscalar density 
we can now write down the relation between M and m in 0(a) improved lattice QCD, 
viz. 

A , M Z A (g )(l + 6 A am q ) 2 

M = ——— -m + Oa , /i=l/L. (6.4 

m(fi) Zp(go,L)(l + bpam q ) 



3 The value of differs from the result quoted in [ p7||33| , because Lmax/n) has been recalculated 
in ref. flSTll . 
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Here, the subtracted mass m q is defined by m q = m,Q — m c , where m c is the critical 
value of the bare quark mass rriQ (c.f. ref. |63||). 

The improvement coefficients 6 a and bp are defined in f33|. A strategy how to 



compute them non-perturbatively has been described recently in 68], but no results for 
our choice of action have been reported so far. However, the difference 6a — 6p, which is 
relevant in the above relation has been shown to be small |3^,[7(]]. Therefore it appears 
safe to neglect 6a and 6p altogether, provided that one is interested in the light quark 
masses only. 

After these considerations we define 

Zm{go) = —-t -Tj—, =z, fi = 1/L, 6.5 

m(/i) Zp(g ,L) 

such that 

M = Z M (go)m(9o) + 0(a 2 ). (6.6) 

Thus the total renormalization Zm consists of the regularization independent part M/fn 
and the ratio Za/Zp, which depends on the details of the lattice formulation. In eq. ( |6.5D 
all reference to the SF scheme has disappeared, so that Zm depends only on the lattice 
regularization. Furthermore Zm, contrary to Zp, is independent of the matching point. 

As described in Sect. [| the ratio M/fn is obtained for scales /i down to a minimum 
value of n = (2L max ) _1 . It is therefore natural to choose 

L = 1.4367-0 (6.7) 



as the value of the matching scale, which is twice the central value of L max in eq. ( |6.1| ). 
We obtain 

M/m(ji) = 1.157(15) at /i = (1.436 rg) ^ . (6.8) 



The error in eq. fl6.8| ) contains a small contribution from the error in eq. ( |6.1| ), which 
was estimated using the perturbative scale dependence of m. 

The next step in the determination of Zm is the calculation of Z^igo) / Zp(go, L) 
at the matching point L = 1.436 ro for a range of bare couplings. Results for the 
normalization constant of the axial current, Z&, have been presented in |6^]. Eq. (6.11) 
in ref. describes Za(<?o) with an accuracy of 1 % for go < 1. 

We have now also computed Zp in the range 6.0 < (3 < 6.5 , (3 = 6/gg. The details 
of the calculation are explained in Appendix |C], and results listed in Table || and in 
eq. (|C.3|). They are represented by a polynomial fit 



Zp(g , L/a) L=1A3ero = 0.5233 - 0.0362 (/3 - 6) + 0.0430 {f3 - 6) 2 , 

f3 = G/gl 6.0</3<6.5, (6.9) 

which is shown together with the data for Zp in Fig. 0. The formula describes the data 
with a statistical uncertainty of about 0.5%. 
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Figure 7: Numerical results for Zp(go, L/a)L=iA36r - The solid line represents the fit eq. ( 3.S ) 
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6.2 
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Table 2: Results for Zp, Za and Zm at typical values of j3 = 6/go 

By combining the results for M/fn, Za and Zp, we can now compute Zm- Results 
at commonly used values of are listed in Table |2|. The entries for Zp and Za in the 
table have been taken from the parametrization eq. Q6.S| ) and Table 2 of ref. ]62fl , respec- 
tively. The first error quoted for Zm comes from the uncertainty in the ratio Za/Zp, 
whereas the second arises from the error in M/rh at the matching scale. It is useful 
to separate the two sources of error, since M/fn has been computed in the continuum 
limit. Therefore, it does not make sense to include the error quoted in eq. (|6.8| ) in the 
continuum extrapolation of lattice data for current quark masses. Instead the uncer- 
tainty of 1.3% in M/rn(fi) at \i = (1.436 ro) -1 should be added in quadrature to the 
quark mass after the extrapolation has been performed. 

As in the case of Zp, it is convenient to represent the results for Zm in terms of 
a polynomial fit function. For bare couplings in the range 6.0 < (3 < 6.5 the data for 
Zyiigo) are well described by the expression 

Z M (g ) = 1.752 + 0.321 (fi - 6) - 0.220 ((3 - 6) 2 . (6.10) 

This parametrization yields Zm with an accuracy of about 1.1%, if the error in M/fn is 
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neglected. 

As explained in detail in Appendix [C] the results for Zp are subject to a systematic 
error of order a, due to imperfect knowledge of the improvement coefficient ct- Our 
estimates of this systematic error show that it is negligible in most applications. In 
fact, the error in the total renormalization Zm is clearly dominated by the uncertainty 
in Za, but not by residual 0(a) effects in Zp. 

We now discuss briefly the application of our results for Zm- Let us assume that 
we want to compute the sum of up and strange quark masses, M u + M s . Using lattice 
data for the hadronic scale tq |3?]] and the bare current quark masses computed in the 
0(a) improved theory one can compute the dimensionless quantity 

Zu{go) x (m u + m s )(g ) x r (g ). (6.11) 

Depending on the chosen values of the bare coupling, Zm can be taken either from Ta- 
ble [2] (using the first of the two quoted errors) or the parametrization in eq. (|6~10l ) (with 
an uncertainty of 1.1%). An extrapolation in (a/ro) 2 to the continuum limit yields the 
value of M u + M s in units of r$. At this point the error of 1.3% in M/fn at the matching 
scale should be added in quadrature to the result. Estimates for the light quark masses 
using our values of Zm will be published elsewhere \71\. 



6.3 Zm for different lattice actions 

The results for Zm presented here have been obtained in 0(a) improved lattice QCD 



using the action defined in [31]. If a different discretization of the QCD action were 
chosen, Zm would have to be recomputed. It is now important to realize that the 
calculation of the universal part of Zm, i-e. the ratio M/fn does not have to be repeated, 
since it is independent of the regularization. This represents a considerable saving of 
CPU time, since the most difficult and demanding part in the determination of Zm is the 
calculation of the scale dependence of M/fn from low to very high energies. Therefore, 
only the calculation of Z^ (as outlined in j62|) and Zp has to be performed for a different 
lattice action. 



7 Conclusions 



In this paper we have demonstrated that the relation between the fundamental pa- 
rameters of QCD and hadronic observables can be computed non-perturbatively with 
controlled errors. Our main results have been obtained in the quenched approximation, 
but our methods carry over literally to full QCD. That is, no additional conceptual 
difficulties are expected, as long as one is interested in QCD with up, down and strange 
quarks only. 

By combining a recursive finite-size technique with lattice simulations we have 
calculated the scale evolution of quark masses in the Schrodinger functional renormal- 
ization scheme, starting from energies of a few hundred MeV up to scales well above 
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100 GeV, where contact with perturbative scaling could be made. In this way we were 
able to extract the ratio M/m at fi = (1.436ro) . Similarly, by computing the scale 
evolution of the running coupling we calculated the A-parameter without compromising 
approximations. By means of a careful removal of lattice artefacts in the computation 
of the scale evolution we have obtained universal results for M/m and the A-parameter. 

Another main result of our work is the calculation of the renormalization factor of 
the pseudoscalar density in 0(a) improved lattice QCD for a range of bare couplings. 
We have computed this factor in the SF scheme at the lowest energy matching point. 
The results can be combined with the previously computed renormalization factor of 
the axial current [S2| and the universal factor M/m. Thereby we have determined the 
total renormalization factor which relates the bare quark masses in the 0(a) improved 
lattice theory to the renormalization group invariant quark masses. 

Our results for A^L and the relation between the bare current quark masses and 
M are not subject to systematic errors other than the use of the quenched approxima- 
tion. In this sense our method to compute renormalization factors of scale-dependent 
quantities is rather unique. Another advantage is that renormalization factors can be 
computed for vanishing quark mass. The method is generally applicable to the problem 
of scale-dependent renormalization and is now also being used to compute renormaliza- 
tion constants for low moments of structure functions [72]. 

Usually quark masses are quoted in the MS scheme at normalization mass \i = 
1 GeV or fi = 2 GeV. Once a non-perturbative solution of the theory becomes possible, 
this convention is not entirely satisfactory, because the MS scheme is only meaningful to 
any finite order of perturbation theory. On the other hand, the renormalization group 
invariant quark masses are non-perturbatively defined and scheme-independent. Hence, 
they are more quotable than the MS masses, and we would like to recommend their use 
in future studies. 



This work is part of the ALPHA collaboration research programme. We are grateful 
to Roberto Petronzio for computer time on the APE / Quadrics at the University of Rome 
"Tor Vergata" , and to Giulia de Divitiis and Tereza Mendes for their help. 
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A Computation of step scaling functions on the lattice 



A.l Calculation of a(u) 

The calculation starts from the lattice step scaling function £(u, a/L) defined as in the 
continuum but for finite resolution a/L. Its precise definition was given in ref. [32]. We 
followed closely this reference. In particular we took exactly the same discretization 



and set the improvement coefficient Ct to its value in one-loop perturbation theory [32], 

c l-l00 P = 1 _ q gg 5 2 _ (A J) 

One then has to be aware that lattice spacing errors linear in the lattice spacing are 
not suppressed completely, because perturbation theory approximates Ct with unknown 
precision. 

Apart from numerical checks and a few simulations on the smallest systems, the 
simulations were performed on APE-100 parallel computers with 128 to 512 nodes. 
We employed the hybrid overrelaxation algorithm described in Sect. 3.1 of setting 
Aqr = L/{2a). For renormalized couplings of around u = 3.5, we used the modified 



sampling introduced in Appendix A of [32] with parameter 7 = 0.02 — 0.05 



Numerical results for running couplings on pairs of lattices L/a, 2L/a at the re- 
spective values of (3 = Q/g^, are listed in Table || for future reference. The data have 
been analyzed as in ref. |j32~[ , by propagating the error of g 2 (L) into a/L), where u 
is the central value of g 2 (L). We then extrapolated a/L) to the continuum limit, 
allowing for the expected dominant discretization error linear in a/L. In practice we fit 

£(u,a/L) = a{u) +uj{u)a/L. (A. 2) 

These fits have good x 2 - We list the fit parameters a(u), the slopes oj(u) and the \ 2 
together with the number of degrees of freedom in the fit, n^f , in Table ||. One observes 
that at the level of our precision lattice artefacts are not significant. Alternatively 
one could assume that ct is well approximated by its one-loop perturbative estimate 
in eq. QA.1| ), so that a/L effects are negligible. An extrapolation using a quadratic 
ansatz, (a/L) 2 , for the leading lattice artefacts would give entirely consistent values for 
cr(u) but smaller statistical errors. 

In Table [l] we have included the results for cr(u) from Table ||] and those in Table 4 



of ref. 1321. 
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Table 3: Pairs of renormalized couplings and the step scaling function S 
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Table 4: Continuum extrapolation of the step scaling function E according to eq. (A. 2) 



A. 2 Calculation of ap(u) 

We now describe the details of the numerical calculations which we have performed to 
determine the step scaling function for Zp. This involves the calculation of correlation 
functions constructed in the framework of the Schrodinger functional in a lattice sim- 
ulation. Here we have used the O(a) improved Wilson action defined in [|63}| . We take 
over the notation introduced in that reference. 



A. 2.1 Definitions 

As explained in Sect. |3], the renormalization constant Zp is defined in terms of correlation 
functions involving the pseudoscalar density and the boundary quark fields. On the 
lattice these correlation functions are represented by 

Mx Q ) = -a 6 ^i^(x) 75 Ir>( 2 ;)C(y)75^ a C(z)), (A.3) 



/l = ~J? £ |(C'(u)7 5 ir a C / (v)C(y)75ir a C(z)). (A.4) 

u,v,y,x 



They have been formally defined in the continuum theory in eqs. ( |3.6| , |3.7| ), and we use 
the same symbols to denote the lattice correlation functions. 

The normalization constant Zp of the pseudoscalar density is then defined by 

Z P (g ,L/a)=c^^ (A.5) 

at vanishing quark mass. The constant c is chosen so that Zp = 1 at tree level. For 
= 1/2, which is chosen in our simulations, we list values for c in Table ||. They differ 
from unity by small terms of order (a/L) 2 . 

We now need to specify the precise definition of the quark mass used in the deter- 
mination of Zp and the step scaling function. The reason is that the point of vanishing 
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current quark mass is not unambiguously defined in a regularization without exact chi- 
ral symmetry. For instance, in the 0(a) improved theory there is an uncertainty of 
order a 2 . The details of the definition will influence the size of 0(a 2 ) lattice artefacts 
in the lattice step scaling function, but not its continuum limit. Following ref. [31| we 
introduce the correlation function of the (unimproved) bare axial current 

/ A (xo) = -a 6 EI^)7075^>(x)C(y)75^ a C(z)), (A.6) 

and define an unrenormalized current quark mass m(go) through 

, s j(go + Wa(sq) + c A ad* dof P (x ) 

m{90) = VA^) • (A - ?) 

Here, do and 8q are the forward and backward lattice derivatives, respectively, and ca 
denotes the coefficient multiplying the 0(a) improvement term in the improved axial 



current. It is understood that eq. (A. 7) is evaluated for 6 = 0, whereas in all other cases 
the correlation functions have been computed for 9 = 1/2. 

The lattice step scaling function Ep is then defined through 

S P (n,a/L) = ^ p( f°' 2L/a) at m(g o ) = 0, g 2 (L)=u. (A.8) 

Here, m(go) is always evaluated for lattice size L/a with xq = L/2. The condition 
m(go) = defines the critical value of the hopping parameter, k c , and the condition 
g 2 (L) = u specifies which value of the bare coupling go is to be used for a given value 
of L/a. 

With these definitions the lattice step scaling function Sp is a function of the 
renormalized coupling u, up to cutoff effects. Its continuum limit is reached as a/L — > 
for fixed u. 

A. 2. 2 Details of the simulation 



The correlation functions in eqs. Q, flP| ) and (|A]|) have been evaluated in a lattice 



simulation using the 0(a) improved Wilson action and the axial current defined in pl[|. 



L/a c L/a c 



6 0.98969662 16 0.99853742 

8 0.99417664 24 0.99934941 

12 0.99740297 32 0.99963393 



Table 5: Values for c, defined implicitly by Zp(0, L/a) = 1 
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In particular, the improvement coefficients c sw and ca given by eqs. (5.16) and (6.5) of 
that reference have been used. 

As in the case of E, boundary 0(a) improvement terms have to be considered. 
In analogy to ct, one also has to cancel the contributions of boundary terms arising 
from local composite operators involving quark fields. This is achieved by including a 
suitably chosen term in the action, whose coefficient is denoted by c~t |33|. Also 2t has 
only been determined in perturbation theory. We have used its perturbative expansion 
to one loop |Tq| 

~l-loop = x _ Qm g g 2 ( A _ 9) 

The influence of the imperfectly known improvement coefficients ct and £t on the contin- 
uum extrapolation has to be investigated. In particular, it is not clear a priori whether 
an extrapolation of Ep in a 2 /L 2 is justified. 



The constants Zp(go, L/a) and Zp(go,2L/a) in eq. (A.8) have been evaluated for 



9 = 1/2 at the point where the current quark mass defined in eq. ( |A.7| ) vanishes. The 
point in k where the condition m(go) = is satisfied is called the "critical" value of the 
hopping parameter, k c . In practice, k c has been obtained through an interpolation of 
the data for m(go) around the point where it vanishes. 

We have evaluated Ep(u, a/L) at each of the 13 chosen values of u for L/a = 6, 8, 12 
and 16. The simulation algorithm is as specified before and the correlation functions 
were computed as detailed in Sect. 2 and Subsect. 3.2 of ref. 1 31 ] . "Measurements" of the 
correlation functions /i , /a and /p were separated by 5 full iterations of the algorithm 
on the smallest lattice, rising to 30 iterations on the largest. We have checked explicitly 
for the statistical independence of our sample by dividing the full ensemble into bins, 
each containing a number of individual "measurements". The statistical errors were 
then monitored as the number of measurements per bin was increased. We did not 
observe any significant change of the errors in Zp for increasing bin size, which we take 
as evidence for the statistical independence of all our samples. 

Since boundary conditions C = C = apply here, the correlation functions /a 
and fp could be averaged with their counterparts f' A and / P , which are defined by a 
time-reflection applied to /a and fp (see eqs. (2.5) and (2.6) of ref. @l). 

The number of measurements was chosen such that the statistical error in Zp(gQ, L/a) 
was typically a factor 2-3 smaller than that of Zp(go,2L/a). Thereby it was ensured 
that the statistical error of Ep was dominated by the uncertainty in Zp(go, 2L/a). Typ- 
ically, we have evaluated Zp(go,L/a) on 200 — 640 configurations, whereas 60 — 120 
measurements were accumulated for Zp(go, 2L/a). 

Statistical errors were computed using the jackknife method. The error in Ep due to 
the uncertainty in the coupling was estimated using the one-loop expansion of Ep in u. 
It turned out to be about 10 times smaller than the statistical error in Ep in the whole 
range of u considered and has therefore been neglected in the error estimate. The error 
in Ep due to the uncertainty in k c has been estimated using the slopes for Zp(go, L/a) 
and Zp(go,2L/a) computed for several values of the hopping parameter. It was found 
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that Ep depends rather weakly on the bare quark mass, so that the uncertainty in Ep 
from the error in k c was negligible. 

Results for Zp(go, L/a), Zp(go,2L/a), k c and the step scaling function Ep are 
shown in Table ||. The data in the last four rows were computed using the two-loop 
formula for ct (see eq. ( |A.11[ ) below). 

A. 2. 3 Continuum extrapolation of Ep 

The next task is the extrapolation of Ep to the continuum limit for fixed coupling u. 
Our results for Ep show at most a weak dependence on the lattice spacing. Before 
performing an extrapolation we first investigate the question whether it is justified 
to assume that the leading cutoff effects are proportional to a 2 /L 2 , given that the 
improvement coefficients Ct and c% are only known in perturbation theory. 

We start by investigating the influence of ct by comparing the step scaling function 
computed using the one-loop expression for 2t to the case where we artificially replace 
the one-loop coefficient by 10 times its value in eq. ( |A.9| ), viz. 

5^ = 1-0.180^. (A. 10) 

For the range of bare couplings used in our simulations, this represents a change of up to 
20%, which is surely a conservative estimate of the remainder of the perturbative series 
for c t . We found that the change in Ep induced by replacing cJ~ 1oo p by c[ amounts to 
around 0.8% for L/a = 6 and u = 3.48, i.e. at the point where the effect is expected 
to be most pronounced. Furthermore, the difference in Ep decreases considerably for 
L/a = 8 and/or smaller couplings. We conclude that the step scaling function is only 
weakly dependent on the value of Ct, and that the effect of using the one- loop estimate 
for ct is negligible at the current level of precision. 

To assess the influence of Ct on the cutoff dependence of Ep, one has to take into 
account that the bare coupling must be adjusted when ct is changed in order to keep the 
renormalized coupling fixed. We have performed this analysis for fixed u = 3.48 where 
the uncertainty in ct is largest. We changed ct by adding the two-loop term, which was 
only known towards the end of our numerical computations [ 35 1 , 



c 2-loop = 1 _ Q8Qg 2 _ o 3 05 4 _ ( A _ H ) 

The corresponding change in Ep is ~ 1.5% for L/a = 6, decreasing considerably as 
the resolution gets finer (cf. Table |6|). We will show in the following that this 0(a)- 
uncertainty, although statistically significant at non-zero a/L, does not affect our values 
for o~p(u) obtained by continuum extrapolations. 
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Table 6: Results for the step scaling function Sp 
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Table 6: (continued) 
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Figure 8: Examples of continuum extrapolations of Ep using Fit B. The dotted lines are the con- 
tinuation of the fit functions to the data points for L/a = 6, which have been excluded from the fit. 



The dependence of Sp on the lattice spacing is modelled using the following fit 
functions 

Fit A : E P (u,a/L) = a P (u) +p{u)a/L (A.12) 

Fit B : S P (u, a/L) = a P (u) + p'{u) a 2 /L 2 . (A.13) 

In other words, we compare the approach to the continuum limit assuming that the 
leading cutoff effects are either linear (Fit A) or quadratic (Fit B) in a/L. Fit A is mo- 
tivated by the fact that not all O(a) improvement terms are known non-perturbatively. 
As a safeguard against higher order cutoff effects, we exclude the data for Sp obtained 
on our coarsest lattices, i.e. for L/a = 6 from the fits. 

The extrapolated values for ap(u) show no significant dependence on the fit ansatz 
(Fit A or B) for all renormalized couplings u. The biggest effect is seen at the largest 
coupling, i.e. u = 3.48, where the results for up (3.48) obtained from either Fit A or B 
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Table 7: Continuum extrapolations of Ep using Fit B 



differ by one standard deviation. Such an effect may be statistical or systematic. To 
test for the latter possibility, we consider the data for £p(3. 48, a/L) obtained using 
c t = c 2 loop . Figure |8| demonstrates that cutoff effects of O(a) are reduced compared to 
Ct = cj~ loop . Furthermore, Fit B applied to the data set for ct = c 2_loop produces an 
extrapolated value which is entirely consistent with the one obtained for ct = Ct _loop . We 
conclude that the small uncertainties which are present in the improvement coefficients 
Ct,Ct are numerically unimportant, and that extrapolations using (a/L) 2 terms as the 
dominant scaling violation are justified. We emphasize that such a statement can only 
be made for a given level of statistical accuracy. Furthermore we use only data for 
L/a = 8, 12 and 16, with the last point already fairly close to the continuum limit. 

For our best estimates we have taken the results from Fit B, obtained for the 
standard one- loop result for ct, and excluding the data for L/a = 6. Typical examples 
of extrapolations at selected values of u are shown in Fig. ||. The fit parameters for all 
extrapolations using Fit B are shown in Table 0, where the point computed using the 
two-loop expression for the improvement coefficient Ct is marked by an asterisk. 



B Error propagation in the scale evolution 

In this appendix we provide further details about the numerical solution of the recursion 
relations (5.3) and ( |5.7| ), the principal aim being to explain how precisely the errors on 
the final results have been determined. 

As discussed in Section B, the calculation starts by fitting the data for the step 
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scaling functions listed in Table |l[ In the case of the function (t(u) our fit ansatz is 

a(u) =u + a u 2 + aiu 3 + . . . + cr n u n+2 , (B.l) 

where 02, . . . , <r n are the fit parameters while 

CT = 21n(2)6 and o\ = a$ + 21n(2)6i (B.2) 

are fixed to the values that one obtains in perturbation theory. The fit function for the 
other step scaling function is taken to be 

ap(u) = 1 + ap t QU + . ■ ■ + (T?,nU n+ , crp,o = — ln(2)d - (B.3) 

Least-squares fits with 1,2 and 3 fit parameters have then been performed and were 
found to represent the data well (the curves shown in Figs. |3] and ||] are for 2 fit param- 
eters) . 

The solution of the recursion relations ( |5.3j ) and ( |5.7| ) is unique once a definite 
expression for the step scaling functions is chosen. In particular, the calculated values 
of Uk and Vk are functions of the fit parameters. However, one should take into account 
that the data for the step scaling functions have errors and thus determine the fit 
parameters only within certain error margins. The errors on the step scaling functions 
are uncorrelated (they come from different simulation runs) and the errors on the fit 
parameters and the u^s and v^s are thus obtained straightforwardly by applying the 
usual rules for error propagation. 

The question now arises to which extent the calculated values and error bounds 
are influenced by the choice of fit functions. We have made two independent checks 
to convince ourselves that this source of systematic error is under control. First we 
noted that compatible results with slightly larger errors are obtained if the number of 
fit parameters is increased. This is the expected behaviour, and since the fit quality is 
already good with one parameter, we decided to quote results for two fit parameters so 
as to be on the safe side. 

The other check that we have made is to apply the fit procedure and error determi- 
nation to simulated data sets that have been generated artificially, assuming cr(u) and 
crp(u) to be equal to some analytically given functions with roughly the right shape. 
Comparing with the "exact" solution of the recursion relations, we then found that the 
systematic bias arising from the choice of fit polynomials is not significant at the current 
level of precision. 

C Computation of Zp at the matching scale 

Here we describe the details of the calculation of Zp(go, £/ a )-L=1.436ro required to con- 
nect the masses in the SF-scheme to the bare current quark masses. Our task is to 
compute Zp for a range of bare couplings, which are commonly used in simulations in 
large physical volumes. In addition, we have to estimate the effect of using perturbation 
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Table 8: Results for Zp(go,L/a) at fixed scale L = 1.436 ro 

theory for the improvement coefficients c t and c t , instead of non-perturbative values, 
which are currently not available. Since ct is known to two loops, we have computed 
Zp using both the one- and two- loop expressions, in order to assess the influence on the 
final estimates. 

Our procedure is as follows. First we choose pairs of ([3, L/a) such that L/a = 
1.436 ro /a, using lattice data for r^/a from |67]]. For a chosen value of L/a this condition 
determines the corresponding (3- value, which is obtained by inserting a/r$ = 1.436(a/L) 
into the parametrization of ln(a/?"o) quoted in eq. (2.18) of ref. 67], viz. 



ln(a/r ) = -1.6805 - 1.7139 (f3 - 6) + 0.8155 {j3 - 6) 2 - 0.6667 (fi - 6) 3 . (C.l) 
Solving numerically for (3 yields the desired combination of L/a). Zp is then com- 



puted for 9 = 0.5 using the normalization condition eq. ( A. 5 ) at the critical value of the 
hopping parameter determined for lattice size L/a. 

This procedure yields (3 with finite accuracy, and the error which propagates into 
Zp can be estimated using one-loop perturbation theory. This error can be neglected, 
since it is more than 10 times smaller than the statistical error in Zp. Results for 
Zp(ga, L/a)L=iA36r computed using both the one- and two-loop expressions for c t are 
listed in Table § 

In principle one could compute Zp for smaller values of (3 and L/a. However, 
for L/a = 6 the condition L = 1.436 ro implies (3 < 6.0, where the parametrizations 



of c sw and ca from ref. |3lJ cannot be used. In order to have an additional value 
at (3 = 6.0, we have applied a special procedure. We have computed Zp for (3 = 
6.0 fixed, but using different lattice sizes for which L/vq straddles the reference value 
L/vq = 1.436. An interpolation in L/tq at fixed (3 then yields the desired result for 

Zp(#o = 1, £/«)£=!. 436 ro- 



Using the parametrization of ro/a, eq. (2.18) in [jpj] we find 

r /a = 5.368(22), = 6.0, (C.2) 
and the results for Zp for three different values of L/tq are listed in Table M. 
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1— loop 
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L/a L/r 


Zp(go = 1, L/a) 


Zp(go = 1, L/a) 


6 1.1177(46) 


0.5728(26) 


0.5712(27) 


8 1.4903(61) 


0.5197(24) 


0.5171(26) 


10 1.8629(76) 


0.4590(31) 


0.4584(30) 



Table 9: Values of Z P at = 6.0 for L around L = 1.436 r 



A quadratic interpolation in L/vq of the data in the table yields 

7 in -1 Tla\ -I °- 5279 ( 24 )' C t" l0 ° P (CV 

Zp(5 °- 1 ' L/a)L=L436r °-\ 0.5253(26), c t 2 - loop ' ^ 

The results for Zp(g , L/ 0)^=1.435 ro f° r 6.0 < < 6.5, using either q _loop or c^ loop 
can be parametrized by a polynomial fit in (/? — 6). For Ct~ loop the result is shown in 
eq. (U). 

We end this appendix with a brief comment on the influence of Ct on the estimates 
of Zp. Similarly to the case of Sp, this has been investigated by repeating the com- 
putation of Zp at (3 = 6.0 using c t with its one-loop coefficient multiplied by a factor 
10. It turned out that the resulting interpolated value of Zp(l, L/a)L=iA36r did not 
change appreciably within statistical errors. We conclude that the influence of ct can 
be neglected at the present level of precision. 
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